3  bulkWES: Tumor Chromosome Deletion

3.1 Set up workspace

# Libraries
library(dplyr)

Attaching package: 'dplyr'

The following objects are masked from 'package:stats':

    filter, lag

The following objects are masked from 'package:base':

    intersect, setdiff, setequal, union

library(tidyverse)

── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ forcats   1.0.0     ✔ readr     2.1.5
✔ ggplot2   3.5.1     ✔ stringr   1.5.1
✔ lubridate 1.9.4     ✔ tibble    3.2.1
✔ purrr     1.0.4     ✔ tidyr     1.3.1

── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag()    masks stats::lag()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

library(ggplot2)

3.2 Load “absolute_seg_file” variable in Terra from all patients

p101_absolute_seg_file <- read.table("/jsimonlab/users/chloetu/melanoma/Melanoma_Eryn/analysis_093024/objects/bulkWES/absolute_seg_file/18279-101.segtab.txt", header = TRUE)
p103_absolute_seg_file <- read.table("/jsimonlab/users/chloetu/melanoma/Melanoma_Eryn/analysis_093024/objects/bulkWES/absolute_seg_file/18279-103.segtab.txt", header = TRUE)
p104_absolute_seg_file <- read.table("/jsimonlab/users/chloetu/melanoma/Melanoma_Eryn/analysis_093024/objects/bulkWES/absolute_seg_file/18279-104.segtab.txt", header = TRUE)
p105_absolute_seg_file <- read.table("/jsimonlab/users/chloetu/melanoma/Melanoma_Eryn/analysis_093024/objects/bulkWES/absolute_seg_file/18279-105.segtab.txt", header = TRUE)
p106_absolute_seg_file <- read.table("/jsimonlab/users/chloetu/melanoma/Melanoma_Eryn/analysis_093024/objects/bulkWES/absolute_seg_file/18279-106.segtab.txt", header = TRUE)
p108_absolute_seg_file <- read.table("/jsimonlab/users/chloetu/melanoma/Melanoma_Eryn/analysis_093024/objects/bulkWES/absolute_seg_file/18279-108.segtab.txt", header = TRUE)
p109_absolute_seg_file <- read.table("/jsimonlab/users/chloetu/melanoma/Melanoma_Eryn/analysis_093024/objects/bulkWES/absolute_seg_file/18279-109_Tumor_WES_01-18279-109_Normal_WES_01.segtab.txt", header = TRUE)
p110_absolute_seg_file <- read.table("/jsimonlab/users/chloetu/melanoma/Melanoma_Eryn/analysis_093024/objects/bulkWES/absolute_seg_file/18279-110_Tumor_WES_01-18279-110-Normal_WES_01.segtab.txt", header = TRUE)
p111_absolute_seg_file <- read.table("/jsimonlab/users/chloetu/melanoma/Melanoma_Eryn/analysis_093024/objects/bulkWES/absolute_seg_file/18279-111_Tumor_WES_01-18279-111-Normal_WES_01.segtab.txt", header = TRUE)

3.3 Format seg files

absolute_seg_file <- do.call(rbind, list(p101_absolute_seg_file, p103_absolute_seg_file, p104_absolute_seg_file, p105_absolute_seg_file, p106_absolute_seg_file, p108_absolute_seg_file, p109_absolute_seg_file, p110_absolute_seg_file, p111_absolute_seg_file))

unique(absolute_seg_file$sample)

[1] "18279-101"                                     
[2] "18279-103"                                     
[3] "18279-104"                                     
[4] "18279-105"                                     
[5] "18279-106"                                     
[6] "18279-108"                                     
[7] "18279-109_Tumor_WES_01-18279-109_Normal_WES_01"
[8] "18279-110_Tumor_WES_01-18279-110-Normal_WES_01"
[9] "18279-111_Tumor_WES_01-18279-111-Normal_WES_01"

# Rename samples and calculate log2ratio
absolute_seg_file <- absolute_seg_file %>%
  mutate(Patient = paste0("P", str_extract(sample, "(?<=18279-).{3}")),
         log2ratio = log2(total_copy_ratio))

# Grab Chr6
chr6_seg <- absolute_seg_file %>%
  filter(Chromosome == 6)

3.4 Plotting LOH for Fig S1C

loh_plot <- ggplot(chr6_seg, aes(x = (Start.bp + End.bp) / 2, y = LOH, color = Chromosome)) +
  geom_rect(fill = "grey", color = NA, alpha = 1, 
            xmin = 29757731, xmax = 33066072, 
            ymin = 0, ymax = 1) +
  geom_segment(aes(x = Start.bp, xend = End.bp, y = LOH, yend = LOH), size = 2, color = "#7b0a0a") +
  labs(title = "CNV Profile Plot of Chr6",
       x = "Genomic Position",
       y = "LOH",
       color = "Chromosome") +
  theme_minimal() +
  facet_wrap(~ Patient) +
  theme(legend.position="none")

Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
ℹ Please use `linewidth` instead.

loh_plot

3.5 Plotting absolute copy numbers of the two alleles for Fig S1D

absolute_cn <- ggplot(chr6_seg, aes(x = (Start.bp + End.bp) / 2, y = modal_total_cn, color = Chromosome)) +
  geom_rect(fill = "grey", color = NA, alpha = 1, 
            xmin = 29757731, xmax = 33066072, 
            ymin = 0, ymax = 7) +
  geom_segment(aes(x = Start.bp, xend = End.bp, y = modal.a1, yend = modal.a1), size = 2, color = "#46786b") +
  geom_segment(aes(x = Start.bp, xend = End.bp, y = modal.a2, yend = modal.a2), size = 2, color = "#a7a486") +
  labs(title = "CNV Profile Plot of Chr6",
       x = "Genomic Position",
       y = "Absolute copy number",
       color = "Chromosome") +
  theme_minimal() +
  facet_wrap(~ Patient) +
  theme(legend.position="none") +
  geom_hline(yintercept = 1, linetype = "dashed", color = "red")

absolute_cn

3.6 Plotting copy number ratio (ie total_copy_ratio column) for Fig S1E

cnr_plot <- ggplot(chr6_seg, aes(x = (Start.bp + End.bp) / 2, y = log2ratio, color = Chromosome)) +
  geom_rect(fill = "grey", color = NA, alpha = 1, 
            xmin = 29757731, xmax = 33066072, 
            ymin = 1.5, ymax = -1.5) +
  geom_segment(aes(x = Start.bp, xend = End.bp, y = log2ratio, yend = log2ratio), size = 2) +
  labs(title = "CNV Profile Plot of Chr6",
       x = "Genomic Position",
       y = "Log2 Copy Number Ratio",
       color = "Chromosome") +
  theme_minimal() +
  facet_wrap(~ Patient) +
  theme(legend.position="none") +
  geom_hline(yintercept = 0, linetype = "dashed", color = "red")

cnr_plot

3.7 Get session info

sessionInfo()

R version 4.3.2 (2023-10-31)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Rocky Linux 8.10 (Green Obsidian)

Matrix products: default
BLAS/LAPACK: /usr/lib64/libopenblasp-r0.3.15.so;  LAPACK version 3.9.0

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
 [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       

time zone: America/New_York
tzcode source: system (glibc)

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] lubridate_1.9.4 forcats_1.0.0   stringr_1.5.1   purrr_1.0.4    
 [5] readr_2.1.5     tidyr_1.3.1     tibble_3.2.1    ggplot2_3.5.1  
 [9] tidyverse_2.0.0 dplyr_1.1.4    

loaded via a namespace (and not attached):
 [1] gtable_0.3.6      jsonlite_1.8.9    compiler_4.3.2    tidyselect_1.2.1 
 [5] scales_1.3.0      fastmap_1.2.0     R6_2.6.1          labeling_0.4.3   
 [9] generics_0.1.3    knitr_1.49        htmlwidgets_1.6.4 munsell_0.5.1    
[13] pillar_1.10.1     tzdb_0.5.0        rlang_1.1.5       stringi_1.8.4    
[17] xfun_0.50         timechange_0.3.0  cli_3.6.3         withr_3.0.2      
[21] magrittr_2.0.3    digest_0.6.37     grid_4.3.2        rstudioapi_0.17.1
[25] hms_1.1.3         lifecycle_1.0.4   vctrs_0.6.5       evaluate_1.0.1   
[29] glue_1.8.0        farver_2.1.2      colorspace_2.1-1  rmarkdown_2.29   
[33] tools_4.3.2       pkgconfig_2.0.3   htmltools_0.5.8.1